Breathing mode for systems of interacting particles. 
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We study the breathing mode in systems of trapped interacting particles. Our approach, based on 
a dynamical ansatz in the first equation of the Bogolyubov-Born-Green-Kirkwood-Yvon (BBGKY) 
hierarchy allows us to tackle at once a wide range of power law interactions and interaction strengths, 
at linear and non linear levels. This both puts in a common framework various results scattered in 
the literature, and by widely generalizing these, emphasizes universal characters of this breathing 
mode. Our findings are supported by direct numerical simulations. 
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Systems of trapped interacting particles are studied 
in many areas of physics: confined plasmas, trapped 
cold atoms, Bose-Einstein condensates, colloidal parti- 
cles, trapped ions, astrophysical systems, the latter ones 
being self confined by the interactions. The low-lying 
oscillatory modes of these systems are a natural object 
of study, as they are an important non destructive tool 
to characterize the system and gain insight into the col- 
lective effects at work. As a consequence, there is an 
abundant literature on the subject, from the different ar- 
eas of research listed above, and corresponding to very 
diverse physical situations: (i) systems with short range 
interactions, such as classical gases or shielded Coulomb 
interaction and (ii) systems with long range interactions, 
such as non neutral plasmas, Coulomb crystals or astro- 
physical system, in which the interactions may be weak 
(gases) or strong (liquids or crystals). 

Diverse approaches and techniques are naturally used 
to investigate these phenomena. A trapped classical gas 
of interacting particles is studied using a Boltzmann- 
Vlasov equation in pQ, where the non linear dynamics 
is approximated with a scaling ansatz, which captures 
the collective effects. Such an ansatz was used earlier for 
the Gross-Pitaevskii equation in [2, 3 . In the confined 
plasma context, the problem is often studied through hy- 
drodynamical equations, in the so-called "cold fluid ap- 
proximation" [4 , where the dispersion relation for fluid 
modes in a cold spheroidal plasma is derived. Follow- 
ing an idea of [5 , Ref. [6 gives an approximate solution 
to the breathing mode of an Id confined plasma beyond 
the cold fluid approximation, using an ad hoc closure 
of the hydrodynamical equations. Monopole modes of 
dusty plasmas interacting with a Yukawa potential are 
investigated in [8]. The breathing mode of trapped 
ions or colloids interacting via Coulomb interactions has 
been studied in Id systems in [9] [10] and in 2d in [IT] 
for crystallized systems, by a direct diagonalization of 
the linearized Newtonian equations of motion. Finally, 
breathing oscillations with attractive interactions have 
been studied in an astrophysical context using the Virial 



theorem [12]. 

Each method applies to a specific situation: Newton 
equations are adapted to a crystallized state with negli- 
gible thermal fluctuations, linearization assumes a small 
amplitude, the Vlasov equation is limited to long range 
interactions and weak correlations. Yet in all cases a sim- 
ilar equation for the breathing mode is obtained. In par- 
ticular, it is intriguing that kinetic descriptions assum- 
ing small correlations between particles, fluid descrip- 
tions, and perturbative expansions around a crystallized 
state all yield similar predictions for the breathing mode, 
at linear and non linear levels. This stunning situation 
calls for a unified theory. In the limit of zero tempera- 
ture, or equivalently infinitely strong interactions, such 
an endeavor has recently been undertaken in the linear 
regime [13]. A more general situation summarizing the 
different possible regimes for a binary isotropic power- 
law interparticle force F(r) ~ l/r k is shown in a diagram 
Fig. [I] We have organized the different cases along two 
axis. On the horizontal axis we represent the interaction 
range, which we will call long-range if k/d < 1 and short 
range otherwise. The case k/d < 1 corresponds to non 
integrable forces at large distances [14 . The vertical axis 
represents the interaction strength with respect to the 
thermal energy. 

In this work, we present a theory of the breathing mode 
of systems of classical trapped interacting particles which 
classifies many cases studied in the references cited above 
in a common framework. The theory is valid both for 
short-range and long-range interactions, for any dimen- 
sion, and for various interaction strengths: it applies to 
the whole diagram of Fig. [I] with restrictions only for 
strongly attractive long range interactions and for even 
moderate attractive short range interactions, where grav- 
itational like collapses in the former case, and strong in- 
stabilities due to the unregularized short range singular- 
ity in the latter are expected. Our theory describes both 
linear as well as non linear oscillations, and isolated sys- 
tems as well as systems in contact with a thermal bath. 

We consider a system of particles confined by an har- 
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FIG. 1: (color online). Diagram of the different regimes for 
the breathing mode. On the horizontal axis, the interaction 
range, measured by k/d, where d is the dimension of the sys- 
tem. The interaction strength is changing along the vertical 
axis. Pictures of some physical examples are inserted for il- 
lustration. 



monic spherical trapping force F trap (r) = — ufir, with 
binary interaction forces ¥i nt . In the canonical setting, 
particles are subjected to a positive constant friction n 
and diffusion D. In the microcanonical setting, n = 0, 
D = and the dynamics is Hamiltonian. To overcome 
the limitations in the validity of the Vlasov equation, 
we describe the cloud of particles by its one-particle 
and two-particles distribution functions /(ri,vi,t) and 
g(i*i, vi, r 2 , v 2 , t). We start from the first equation of the 
BBGKY hierarchy, which we complement by a Fokker- 
Planck operator to include the temperature in the canon- 
ical case: 



of 

dt 



+ V r .(v/) + F trap .V v / + C[g] = DA v f + «V v .(v/), 

(1) 



where C[g] is the interaction term given by: 

C[#](ri, vi,t) = J F in£ (ri,r).V Vl ^(ri, vi,r, v,t) drdv. 

(2) 

We stress that Eq. ([!]), in contrast with the Vlasov equa- 
tion, can also describe strongly correlated systems. We 
assume in the following the existence of a stationary state 
fo and go, not necessarily the thermodynamic equilib- 
rium [15] . We now drastically simplify the dynamics by 
using a scaling ansatz (TJ [2j [3] , which we extend here to 
the two-particles function g: 



/(ri,vi,t) = /o00i,vi)) 
#0i,vi,r 2 ,v 2 ,t) = 0oW>(ri,vi,r 2 ,v 2 )) 



(3) 



parameter A. Introducing Eq. (|3| into Eq. leads to: 
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-DA 2 ^} + C[g o^](n, vi,t) = 0, 



(4) 

where the difficulty is to deal with the interaction term. 
We now assume a homogeneous two-body interaction 
with degree —k: 



F int (Ari,Ar 2 ) = ^F int (ri,r 2 ), 



(5) 



as for example a pure power law. The important step is to 
replace the interaction term C[gooiJ;](ri : vi, t) by a linear 
combination of fo and its derivatives. This is achieved 
using the condition ^ and the fact that fo and go are 
stationary solutions of Eq. ([!]). Equation Q becomes 
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Multiplying the previous equation by RjVj/N, and in- 
tegrating over dUdY, we obtain a constraint on the pa- 
rameter A: 



A + kX + A 
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(R 2 j)fo 



= 0, (7) 



where j is a coordinate label, and we have set (x)f = 
J ^(r, v)/(r, v, t)drdv. In the dynamical equation for 
A (Eq. ([7])), all parameters are computed as averages over 
the stationary distribution fo . For Eq. ^ to be a unique 
equation, it is necessary that the ratio (V 2 ) f / (R 2 ) / 
does not depend on j, which is true if the trap and in- 
teractions are isotropic. 

We introduce the dimensionless parameter p = 
(V?)f /(w$(Rj)fo) ~ k B T/E trap , wher k B T is the ther- 
mal energy and E trap the typical potential energy due to 
the trap. At the canonical equilibrium, (V 2 )f = uJqL 2 , 
where L is the typical size of the system without inter- 
action. The parameter p = L 2 / (R'j) / thus describes 
change of the square of the size of the trap due to the in- 
teractions. The range p < 1 (resp. p > 1) corresponds to 
a repulsive (resp. attractive) interaction. A value of the 
parameter p ~ 1 means high temperature or negligible 
interactions. The limits p — > and p — > +00 correspond 
to zero temperature or strong repulsive and attractive 
interaction. We can now rewrite Eq. ^ as 



with y?(ri,vi) = (Ri = ri/A, Vi = Avi — Ari) and 
^(ri, vi,r 2 , v 2 ) = (<p(ri, vi),<p(r 2 , v 2 )). All time depen- 
dence in the dynamics is now included in the positive 



A + «A + 0'(A) =0, (8) 
which corresponds to the equation of a damped anhar- 
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monic oscillator in the potential <f>: 

MX) = { ^ (** + 5* + ^ Al_fe ) ' iffc ^ X ' 
\ Wo 2 (iA 2 + if + (p-l)logA), if fc = 1 . 

(9) 

The first term in Eq. ([9| is the quadratic confining poten- 
tial, the second one corresponds to a pressure term and 
the last one is introduced by the two-body interaction. 

For repulsive interactions (p < 1), <\> is strictly convex 
for all k ^ 0. It diverges as A -2 when A — > and as 
cJqA 2 /2 when A — > +oo. Its unique minimum is A = 1. 
The A -2 divergence at small A is due to pressure effects 
for very compressed clouds, and thus does not depend on 
the interaction. It yields a generic shape for the breath- 
ing oscillations in the nonlinear regime. For attractive 
interactions (p > 1), if < k ^ 3, <fi has exactly the same 
qualitative properties as in the repulsive case. For k > 3, 
(j) tends to minus infinity when A goes to zero, indicating 
a possible collapse of the cloud. However, such unregu- 
larized power law forces may create instabilities, so that 
these predictions for the breathing dynamics should be 
considered with care. 
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FIG. 2: (color online). Frequency of the linearized breathing 
mode as a function of the interaction strength p, for different 
values of interaction range k. 

From Eq. ([9]), we obtain the general expression of the 
breathing oscillation frequency in the small friction limit, 
as a function of the interaction range k and the interac- 
tion strength p: 

u;(k,p)=u;oi(3-k)(p-l)+4} 1/2 . (10) 

This expression recovers the well known limits uo = 2o;o 
for a non interacting gas (p = 1) and uo = V3ujq for a 
strongly interacting Coulomb plasma (p = 0, k = 2). It 
provides a generalization to the whole (k,p) plane shown 
in Fig. [2] and is independent of the dimension. We note 
that in 3 dimensions, the breathing frequency is a de- 
creasing function of the interaction strength for repulsive 



long range interactions, and a increasing function of the 
interaction strength for repulsive short range interaction. 

We can now compare the general Eq. ^ to the results 
found in the literature for various specific situations. Os- 
cillations of crystallized systems [9j [TU1 HTJ [13] correspond 
to negligible pressure effects, i.e. p — and the A -3 term 
of is absent. In [6 , the authors consider a Id plasma 
(k = 0) with p not too small, and introduce a pressure 
yielding the A -3 term, which leads to the exact equivalent 
of Eq. ([7]) . Note that Eq. ([7]) also contains the case of a 
classical gas with "mean field" interactions pQ. This work 
considers a Dirac 5 potential, which corresponds to a ho- 
mogeneity degree —k = —d— 1. This result emphasizes 
that the present theory is not only valid for power-law 
forces. 

In order to test the domain of validity of the ansatz so- 
lution, we have performed numerical simulations for dif- 
ferent force index /c, parameter p and amplitude of initial 
perturbation, in two and three dimensions, with (canon- 
ical ensemble) or without (microcanonical ensemble) a 
thermostat. We simulate the system using a molecular 
dynamics approach with N = 4000 particles. The inte- 
grator scheme is a Verlet-leapfrog algorithm [16] in the 
microcanonical or canonical ensemble. The forces are ex- 
actly computed at each time-step. As strong short range 
singularities for parameters in the upper right corner of 
Fig. [T] create numerical difficulties, we have not tested 
the theory in this region. The computer simulations are 
performed as follows: we first equilibrate the system in 
a stationary state fo. Then, at t = 0, we introduce a 
perturbation by rescaling the positions and velocities ac- 
cording to Eq. |3| and we let the system evolve. A similar 
simulation of a Id Coulomb system in the microcanon- 
ical ensemble has been performed in [6]. The results of 
our extensive simulations may be summarized as follows: 
(i) Eq. ([7| always picks up quite precisely the oscillation 
frequency, but not always the amplitude decay, (ii) For 
strongly repulsive interaction (p — > 0), Eq. ^ describes 
very precisely the whole dynamics. (Hi) For a repulsive 
long-range or short range interaction and intermediate 
p (i.e. p ~ 0.5) the agreement for the oscillation am- 
plitude is not perfect (Fig. [3|. (iv) For attractive long 
range interactions, the accuracy of the ansatz degrades 
as p increases (Fig. [4|. 

To explain these results, we first stress that in the limit 
p — > 0, Eq. ([7]) is exact. In this case, it may indeed be de- 
rived directly from Newton equations, as done in [13] in 
the linear approximation. The correct generalization for 
an arbitrary perturbation amplitude is given by Eq. 
For intermediate p, we attribute the discrepancy between 
the predicted and simulated oscillation amplitudes to ef- 
fects that are not taken into account in the simple dy- 
namical ansatz Q , and thus limit the validity of Eq. ^ . 
Indeed, for long range interactions, one would expect col- 
lective effects neglected in the ansatz (Landau damping, 
phase mixing, etc) to play a role in the oscillation decay 
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beyond the friction k). Similarly, for short range inter- 
actions, two-body collisions may be important. This ex- 
planation is supported by the frictionless microcanonical 
simulations: when there is no amplitude decay in the mi- 
crocanonical ensemble, which means that phase mixing 
and two body collisions are negligible, Eq. ([7]) correctly 
predicts the breathing frequency and amplitude, with or 
without friction. Conversely, amplitude decay or modu- 
lation in the microcanonical ensemble is associated with 
discrepancies between theory and simulations. 
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FIG. 3: Evolution of the typical size of the cloud, in a case 
where the ansatz does not describe the full dynamics. The 
space dimension is d = 2, and the interactions are repulsive. 
The parameters are k = 4 (short range interaction), ooo/k — 
17.8 and p = 0.63. 



BBGKY hierarchy and a scaling ansatz for the dynam- 
ics, we have derived a non-linear equation describing 
the breathing oscillations of trapped particles interact- 
ing via homogeneous forces. The derivation and equa- 
tion are valid independently of the temperature, inter- 
action strength, interaction range and dimensionality of 
the physical space, and it is successfully confronted to di- 
rect numerical simulations. The main limitation is due to 
phase mixing phenomena for long range interacting sys- 
tems and two-body collisions in short range interacting 
ones, especially for weak repulsive and attractive inter- 
actions, where they introduce damping and loss of co- 
herence, unaccounted for in the scaling ansatz. Even 
though we have concentrated on power-law interactions, 
the homogeneity condition for the force is more general. 
It includes for instance dipolar interactions and some non 
potential forces, allowing the use of the ansatz technique 
in such cases. Finally, the general results obtained in 
this letter might be a useful guide for experimentalists to 
extract information on the interaction between particles 
from easily measurable phenomena, such as breathing 
mode dynamics. 

Beyond the breathing mode, a generic study of 
quadrupolar modes would be very desirable, as harmonic 
traps are often anisotropic in experimental situations. 
This is not possible with the scaling ansatz, except in 
special cases. However, following the lines of the present 
letter, and applying methods used in [8], a more general 
approach should be possible. 




25 50 

UJ t 



FIG. 4: Evolution of the typical size of the cloud, in a case 
where the ansatz does not describe the full dynamics. The 
frequency is still very well captured. The space dimension is 
d — 3, and the interactions are attractive. The parameters 
are k = (long range interaction), ujo = 17.8, k = (micro- 
canonical ensemble) and p = 2.2. Same parameters for the 
inset except p = 70. 

In summary, starting from the first equation of the 
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